# Replication code for Tables and Figures for James Fearon, "Cooperation, Conflict, and the Costs of Anarchy."

# First prepare country year level arms data.
setwd("../papers/coopconf/IO_final")
d = read.csv("cyear_raw.csv",as.is=T)
d = d[order(d$ccode,d$year),]

 # define measures of military effort.
 # note:  gdp.currentUS is from World Bank, WDI 2015.  I have added population estimates for Taiwan,
 # (which is not a WB member) taken from the IMF World Economic Outlook.

 # Inspection also revealed some minor errors in the COW milex series in NMC 4.0.
 # Vietnam 2007 appears to be data input error.  I multiply this by 1000 to get (presumably) correct value.
 # Cuba 1994 also appears to be a mistake. Set to NA.

  # look at milex series to correct a few errors.
d$milex[d$ccode==816 & d$year==2007] = 3709000  # looks like data entry error.
d$milex[d$ccode==40 & d$year == 1994] = NA  # Cuba.

 # compute COW based military burden.
d$milbur = d$milex*1000/d$gdp.currentUS
 
 # soldiers per 1000.
d$milpc = d$milper*1000000/d$pop

 # Use SIPRI military burden estimates (which are very similar to COW's), to add data for 2008-2010.
d$milbur[d$year > 2007] = d$sipri.milbur[d$year > 2007]
  # fill in a few other missing COW milex years with SIPRI.  Have checked that this is good.
i = d$year > 1987 & d$year < 2008
table(is.na(d$milbur[i]),is.na(d$sipri.milbur[i]))  
d[i & is.na(d$milbur) & !is.na(d$sipri.milbur),c('country','ccode','year','milbur','sipri.milbur')]  
d$milbur[i & is.na(d$milbur) & !is.na(d$sipri.milbur)] = d$sipri.milbur[i & is.na(d$milbur) & !is.na(d$sipri.milbur)]

# Table 1.

 # Bring in Thompson and Dreyer rivals data.
rival = read.csv('Rivalry data 1816-2010 with type.csv',as.is=T)
 # prepare data
rival$start[rival$start=="pre1816"] = 1816
rival$end[rival$end=="ongoing"] = NA
rival$start = as.numeric(rival$start)
rival$end = as.numeric(rival$end)
rival[,7:10] = apply(rival[,7:10],2,function(x) replace(x,is.na(x),0))
rival[order(rival$ccode1,rival$ccode2,rival$start),]
 # consider post 1945 for this Table.
rival = rival[rival$end > 1945|is.na(rival$end),]

 # directly count types of rivalry for each state in these data for Cold War and post-Cold War.
rival$end2 = rival$end
rival$end2[is.na(rival$end2)] = 2010  # note T and D data end in 2010 so use only d$year < 2011 below
 
c = sort(unique(c(rival$ccode1,rival$ccode2)))
crival = data.frame(ccode = c)
crival$terr.CW = sapply(c,function(x) nrow(rival[rival$Spatial==1 & (rival$ccode1==x|rival$ccode2==x) & rival$start < 1991,]))
crival$terr.pCW = sapply(c,function(x) nrow(rival[rival$Spatial==1 & (rival$ccode1==x|rival$ccode2==x) & rival$end2 > 1990,]))
crival$other.CW = sapply(c,function(x) {i = rival$Spatial==0 & (rival$ccode1==x|rival$ccode2==x) & rival$start < 1991
	nrow(rival[i,]) # note we are counting number of states with which c has a non-spatial rivalry
	})
crival$other.pCW = sapply(c,function(x) {i = rival$Spatial==0 & (rival$ccode1==x|rival$ccode2==x) & rival$end2 > 1990
	nrow(rival[i,]) # note we are counting number of states with which c has a non-spatial rivalry
	}) 

 # Cold War and post Cold War.  top military burdens by decade, using states with population over .5 million
 # also compute milpc at same time.
m = 499999 
a = d$ccode[d$year==1990 & d$pop > m]  # but what biggish states entered after 1990?
b = d$ccode[d$yrentry == d$year & d$yrentry > 1990 & d$pop > m]
z = d[d$ccode %in% c(a,b) & d$year < 2011,]
a = data.frame(tapply(z$milbur,list(z$ccode,z$year > 1990),mean,na.rm=T))
t  = tapply(z$milpc,list(z$ccode,z$year > 1990),mean,na.rm=T)
a$milpc.CW = t[,1]; a$milpc.pCW = t[,2]

a$country = d$country[match(rownames(a),d$ccode)]
a$ccode = as.numeric(rownames(a))
a = cbind(a,crival[match(a$ccode,crival$ccode),2:5])
a[,7:10][is.na(a[,7:10])] = 0
names(a)[1:2] = c('mb.CW','mb.pCW')

 # milbur part of table.  
z = paste(a$country[order(-a$mb.CW)][1:15],"$^{",a$terr.CW[order(-a$mb.CW)][1:15],",",a$other.CW[order(-a$mb.CW)][1:15],"}$",sep="")
z2 = paste(a$country[order(-a$mb.pCW)][1:15],"$^{",a$terr.pCW[order(-a$mb.pCW)][1:15],",",a$other.pCW[order(-a$mb.pCW)][1:15],"}$",sep="")
z3 = paste(a$country[order(a$mb.CW)][1:15],"$^{",a$terr.CW[order(a$mb.CW)][1:15],",",a$other.CW[order(a$mb.CW)][1:15],"}$",sep="")
z4 = paste(a$country[order(a$mb.pCW)][1:15],"$^{",a$terr.pCW[order(a$mb.pCW)][1:15],",",a$other.pCW[order(a$mb.pCW)][1:15],"}$",sep="")

z5 = rbind(cbind(z,z2),cbind(z3,z4))
 
 # milpc part of table.  
z6 = paste(a$country[order(-a$milpc.CW)][1:15],"$^{",a$terr.CW[order(-a$milpc.CW)][1:15],",",a$other.CW[order(-a$milpc.CW)][1:15],"}$",sep="")
z7 = paste(a$country[order(-a$milpc.pCW)][1:15],"$^{",a$terr.pCW[order(-a$milpc.pCW)][1:15],",",a$other.pCW[order(-a$milpc.pCW)][1:15],"}$",sep="")
z8 = paste(a$country[order(a$milpc.CW)][1:15],"$^{",a$terr.CW[order(a$milpc.CW)][1:15],",",a$other.CW[order(a$milpc.CW)][1:15],"}$",sep="")
z9 = paste(a$country[order(a$milpc.pCW)][1:15],"$^{",a$terr.pCW[order(a$milpc.pCW)][1:15],",",a$other.pCW[order(a$milpc.pCW)][1:15],"}$",sep="")

z10 = rbind(cbind(z6,z7),cbind(z8,z9))

z11 = cbind(z5,z10)

 # TABLE 1
apply(z11,1,function(x) cat(cat(x,sep=" & "),"\\\\\n"))


 # Figure 2 code, territorial rivalries over time.

 #RIVALS

#opening datasets

rival = read.csv('Rivalry data 1816-2010 with type.csv',as.is=T)
nmc <- read.csv("NMC_v4_0.csv",as.is=T,na="-9")

#setting up rival data
rival$start <- as.character(rival$start)
rival$end <- as.character(rival$end)
rival$start[rival$start=="pre1816"] = "1816"
rival$end[rival$end=="ongoing"] = "2015"
rival$start = as.numeric(rival$start)
rival$end = as.numeric(rival$end)
rival[,7:10] = apply(rival[,7:10],2,function(x) replace(x,is.na(x),0))

# making dataset by ccode for all rivalries.

r = rival
r2 = r
names(r2) = c('ccode2','abbrev2','ccode1','abbrev1','start','end','Positional','Spatial','Ideol','Interv')
r2 = r2[,c(3,4,1,2,5:10)]
r = rbind(r,r2)
r = r[order(r$ccode1,r$ccode2,r$start),]

# making dataframe for 1816:2010

a = data.frame(year=1816:2010)
a$terr = sapply(a$year,function(x) sum(x >= r$start[r$Spatial==1] & x <= r$end[r$Spatial==1]))
a$oth = sapply(a$year,function(x) sum(x >= r$start[r$Spatial==0] & x <= r$end[r$Spatial==0]))

# counting rivalries per state based on the COW list
nmc$rival.terr = 0  # number of territorial rivalries for country year
for (i in 1:nrow(r)) {j = nmc$ccode==r$ccode1[i] & nmc$year >= r$start[i] & nmc$year < r$end[i]
nmc$rival.terr[j] = nmc$rival.terr[j] + r$Spatial[i]	
}

t = tapply(nmc$rival.terr,nmc$year,mean,na.rm=T) 
a$pctterr = t[match(a$year,names(t))]  # mean number rivalries per state by year.

# counting rivalries for the 15 % most powerful states, using COW CINC score for "power".
t = tapply(nmc$cinc,nmc$year,function(x) quantile(x,.85,na.rm=T))   
nmc$cinc85 = t[match(nmc$year,names(t))]
i = nmc$cinc > nmc$cinc85
t = tapply(nmc$rival.terr[i],nmc$year[i],mean,na.rm=T)
a$pctterr.top15pct = t[match(a$year,names(t))]

pdf("terrrivals.pdf") 
plot(pctterr.top15pct ~ year,data=a,type='n',col=1,ylim=c(0,4),xlab="year",ylab="territorial rivalries per state")
rect(1914,-1,1918,5,col="gray",lty=0)
rect(1939,-1,1945,5,col="gray",lty=0)
lines(pctterr.top15pct ~ year,data=a,col="black",lwd=2)
lines(pctterr ~ year,data=a,col=1,lty=2,lwd=2)
legend(1820,.6,c('all states','top 15% on military capabs \n (COW CINC scores)'),lty=2:1,lwd=2,col=1,cex=.9)
#title("Figure 2. Territorial rivalries per state, 1816-2007")
graphics.off() 
 
 
# FIGURE 3 (major power arms burdens over time) 
 
 # bring in COW wars and compute 
cow = read.csv('InterstateWarData_v4.0.csv',as.is=T) 

 # add war markers to nmc
nmc$iwar = 0
for (i in 1:nrow(cow)) {j = nmc$ccode==cow$ccode[i] & nmc$year >= cow$StartYear1[i] & nmc$year <= cow$EndYear1[i]
	nmc$iwar = nmc$iwar + j
	}
a = cow[cow$StartYear2 > 0,]
for (i in 1:nrow(a)) {j = nmc$ccode==a$ccode[i] & nmc$year >= a$StartYear2[i] & nmc$year <= a$EndYear2[i]
	nmc$iwar = nmc$iwar + j
	}

maj = c(2,200,220,255,325,365,640,710,740)   # Major powers in Figure.  Omitting AUH because no post 1918.
nmc$ccode[nmc$ccode==260] = 255  # assign West Germany to Germany post WWII.

 # DV and omit war years.
nmc$milpc = 1000*nmc$milper/nmc$tpop  # Soldiers per 1000.

nmc$milpc2 = nmc$milpc  
nmc$milpc2[nmc$iwar>0] = NA  # this will be the solid line, peacetime spending.
nmc$milpc2[nmc$year %in% c(1918:1921,1945:1947)] = NA # dropping immediate demobilization years of the WW's.
nmc$milpc2[nmc$ccode==2 & nmc$year %in% c(1860:1865)] = NA  # US Civil War, as noted in text
nmc$milpc3 = nmc$milpc  
nmc$milpc3[nmc$iwars ==0] = NA  # for dotted line indicating war years.

library(locfit)  # smoother
nmc$statenme[nmc$ccode==2] = "U.S.A."

pdf("milpc_major_powers.pdf",width=9) 
par(mfrow=c(3,3),oma = c(0,0,1,0),mar = c(3,4,2,1)) #,mai = c()) 
for (i in maj) {j = nmc$milpc2[nmc$ccode==i][match(1816:2007,nmc$year[nmc$ccode==i])]
				plot(1816:2007,j,'l',xlab="",ylab="armed forces/1000")
				title(unique(nmc$statenme[nmc$ccode==i][1]),line=.5)
				title(xlab = "year",line = 2)
				lines(1816:2007,nmc$milpc3[nmc$ccode==i][match(1816:2007,nmc$year[nmc$ccode==i])],lty=2)
				lines(locfit(j ~ lp(I(1816:2007),nn=.5)),lty=3,lwd=1.5,col=1) }				
#title("Figure 3.  Peacetime armed forces/1000 population",outer=T)
graphics.off()


